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Abstract 

Controlling the spatial distribution of liquid droplets on surfaces via surface energy patterning 
can be used to control material delivery to specified regions via selective liquid/solid wetting. While 
studies of the equilibrium shape of liquid droplets on heterogenous substrates exist, much less is 
known about the corresponding wetting kinetics. We make significant progress towards elucidating 
details of this topic by studying, via large-scale atomistic simulations, liquid nanodroplets spreading 
on chemically patterned surfaces. A model is presented for lines of polymer liquid (droplets) on 
substrates consisting of alternating strips of wetting (equilibrium contact angle 0q ~ 0°) and non- 
wetting (6q ~ 90°) material. Droplet spreading is compared for different wavelength A of the 
pattern and strength of surface interaction on the wetting strips. For small A, droplets partially 
spread on both the wetting and non-wetting regions of the substrate to attain a finite contact 
angle less than 90°. In this case, the extent of spreading depends on the interaction strength in 
the wetting regions. A transition is observed such that, for large A, the droplet spreads only on 
the wetting region of the substrate by pulling material from non-wetting regions. In most cases, a 
precursor film spreads on the wetting portion of the substrate at a rate strongly dependent upon 
interaction strength. 
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I. INTRODUCTION 



The spreading of a liquid droplet on a surface proceeds as the droplet balances the 
interfacial tensions between the solid, liquid, and vapor phases. On a homogenous substrate, 
completely wetting droplets (i.e. where the equilibrium contact angle #0 = 0°) often form 
a monolayer film or precursor foot on the surface that spreads ahead of the main droplet. 
This has been observed in experiments as well as simulation DDDDO . This 

foot grows diffusively with an effective diffusion constant D e R that depends on the square 
root of the initial radius R Q of the droplet Q|. While much has been revealed about the 
kinetics of droplets spreading on homogenous surfaces, spreading on chemically patterned 
surfaces is much less understood. In particular the influence of the pattern feature size on 
the spreading rate of both the precursor foot and main droplet is unknown. 

With the advent of modern lithography techniques, it is now possible to chemically pat- 
tern a substrate with lyophobic and lyophilic surface domains on the nano to micron scale 
. These chemically structured surfaces can be useful in controlling the wetting of liquid 
droplets. Depending on the size of the pattern and the contact angle in the two domains, 
wetting transitions can occur in which the liquid droplet can change its shape or morphol- 
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The equilibrium shapes of liquid droplets on a variety of chemica 
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surfaces has been studied both experimentally and theoretically 
Introduction of a lyophobic grid can be used to confine liquid droplets and has been shown 
to improve the quality of ink jet printing [21]. In most of these studies, the liquid droplet is 
confined by the pattern of lyophilic and lyophobic domains and does not spread. 

An alternate situation occurs if one places a long liquid line transverse to a chemically 
striped surface of alternating lyophilic and lyophobic domains. This is the situation we study 
here for the case of liquid polymer droplets of varying chain length. For small A, the main 
droplet partially wets both regions of the surface equivalently and does not spread even in 
the wetting region. For sufficiently large A the main droplet spreads on the wetting region of 
the substrate by pulling material from non-wetting regions. In most cases, a precursor film 
spreads on the wetting portion of the substrate at a rate strongly dependent upon interaction 
strength even when the main droplet does not spread. 
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II. COMPUTATIONAL PROCEDURE 



Here we present extensive molecular dynamics (MD) simulations of the spreading of 
short chain polymer droplets spreading on a chemically patterned surface. The polymers 
are modeled by a coarse grained bead spring model in which each polymer contains N beads 
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per chain with N = 10 and 100. As the entanglement length N e ~ 72 for this model 
chains are unentangled. The polymer is represented by spherical beads of mass m attached 
by springs, which interact with a truncated Lennard- Jones (LJ) potential, 
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where e and a are the LJ units of energy and length and the cutoff r c = 2.5a. The monomer- 
monomer interaction e is used as the reference and all monomers have the same diameter 
a. For bonded monomers, we apply an additional potential where each bond is described 
by the FENE potential |23J with k = 30 e/a 2 and Rq = 1.5 a. The substrate is modeled as 



a flat surface since it was found previously 



that with the proper choice of thermostat, 



simulations using a flat surface exhibit the same behavior as a realistic atomic substrate 
with greater computational efficiency 

The interactions between the surface and the monomers in the droplet at a distance z 
from the surface are modeled using an integrated LJ potential, 
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where e w is the monomer- surface interaction strength and z c = 2.2a. Extending the range of 
this surface interaction to infinity shifts the wetting transition to a lower energy; however, 
the qualitative spreading behavior is identical to the z c = 2.2a case ^(J. For z c = 2.2a, 
the critical wetting strength (i.e. 9 — > 0°) is e c w ~ 1.75s for iV = 10 and e c w ~ 2.25e for 
iV = 100 [9|. The system is periodic in the y direction with length L y and open in the other 
two directions. The chemical pattern consists of infinite strips of wetting and non-wetting 
regions of equal width A/2 as shown in Figs. 1-2 and Fig. 4. In the present simulations the 
monomer/wall interaction e w depends on the position of each monomer in the y direction. 
For the non-wetting region, e„ was chosen so that 8q ~ 90°, which for N = 10 is e w = 0.9e 
and for iV = 100 is e w = l.Oe |9(. In the wetting region, two cases were studied for iV = 10: 



e w = 2.0e and e w = 3.0e, corresponding to weak and strong wetting, respectively. For 
N = 100, e w = 2. be, which is in the weak wetting regime. All of the droplets presented here 
are modeled as hemicylinders as described previously [9] with initial droplet radii R = 50a. 
The effect of varying Rq for homogeneous droplets has been studied previously [10]. 

To provide a realistic representation of the transfer of energy in the polymer droplet, a 
Langevin thermostat is applied only near the surface, 

rrii— i = -AUi - miL-^ + W<(*), (3) 

where to, is the mass of monomer i, is the position of monomer i, 7x is the friction 
coupling term for the Langevin thermostat, —AUi is the force acting on monomer i due to 
the potentials defined above, and Wj(t) is a Gaussian white noise term. We use a Langevin 
coupling term with a damping rate that decreases exponentially away from the substrate 



24( | . 7l(z) = Yl ex P ( a ~ z )- Here Yl = 3.0r _1 and z is the distance from the substrate. The 
effect of varying 7£ for homogeneous surfaces has been discussed in ref. [8j]. By including 
only a coupling of monomers near the surface to the thermostat, we avoid the unphysical 
effect of screening the hydrodynamic interactions in the droplet. We also damp monomers 
near the surface stronger than those in the bulk, which is a more physically reasonable 
description of dissipation mechanisms due to the surface interaction. 

In the simulations L y varied from 60a to 200<j, depending on the periodicity of the pattern 
A. For A < 60er, L y = 60cr while for A > 100a, L y = A. Results for a homogenous wetting 
substrate, labeled A = oo are for L y = 60a. The total size of the droplets varied from 
iV = 200 000 to 680 000 monomers depending on the value of L y . The vapor pressure is low 
so that spreading does not occur via vaporization and condensation. However for small A, 
particularly for iV = 10, a few chains separate from the precursor foot and randomly diffuse 
on to the non-wetting region of the substrate and sometimes evaporate from the substrate. 
These chains are subsequently removed from the simulation once they move away from the 
surface a distance several times Rq. 

The equations of motion are integrated using a velo city- Ver let algorithm with a time 

. The simulations are performed at a temperature 
, modified to include the chemically patterned surface. 
The instantaneous contact radius for the main droplet r b (t) and the radius of the precursor 
foot rf(t) are determined every 400r according to the procedure described previously 0. 



Briefly, the system is divided into layers of thickness 1.5<r parallel to the interface. Data 
presented for foot kinetics are obtained from the layer closest to the surface (layer 1) and 
data for bulk kinetics are obtained from layer 4, which is far enough above the substrate to 
provide an accurate description of bulk droplet behavior. 

III. RESULTS 

Time sequences for N = 10 and iV = 100 are shown in Fig. 1 for A = 200a. For the 
main droplet to spread (bulk spreading), the gain in interface energy must exceed the cost 
of increasing the surface area of the droplet. This balance depends on A, N, and e w . For 
N — 10, bulk spreading occurs on the lyophilic strip as material is transfered from the non- 
wetting to the wetting region. As such, the cross-section of material in the non-wetting region 
shrinks while maintaining a hemi-cylinder shape with contact angle near 90°. At later time in 
the simulation, material is completely removed from non-wetting regions. Given periodicity 
in the y direction, this demonstrates a case where an intial line of liquid (in y) transforms 
into a series of parallel liquid strips (in x) with separation distance in y equal to A/2. For 
N = 100, evidence of N dependence is observed in that bulk spreading does not occur on 
the wetting region. A precursor film or foot is observed to spread ahead of the main droplet 
on the lyophilic strip in both systems. Thus, even when the main droplet does not spread, a 
precursor foot may spread on the wetting region of the substrate. For all cases where a foot 
was observed to spread, the contact radius follows rj(t) = 2D cS t where D eS is the effective 
diffusion coefficient. For N = 10 and A = 200a (Fig. 1), the effective diffusion constant 
D e g of the precursor foot is essentially equal to that for the pure wetting case A = oo. 
Spreading of the precursor foot for N = 100 causes depletion of material on the wetting 
region. Furthermore, bulk material is unable to move quickly enough from the non-wetting 
region to the wetting region, resulting in the formation of a saddle shape at the center of 
the droplet. At late times, the precursor film spreads by polymer moving perpendicularly 
from the non-wetting to the wetting region, then continuing along the wetting region. In 
this case D cS is 40% of that for A = oo. Although we do not access such time scales in our 
simulations, we expect that material will completely wick away from the non-wetting region 
as the precursor film continues to spread. 

Figure 2 presents snapshots of iV = 10 droplets at late time for varying A. The snapshots 
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FIG. 1: (color) Sequence of simulation snapshots for the N = 100 (left) and N = 10 (right) 
droplets for A = 200a for (a) t = 0,(b) t = 12, 000, (c) t = 24, 000, (d) t = 36, 000, (e) t = 48, 000r. 
For N = 100, in non-wetting regions e w = l.Oe (green), while in wetting regions e w = 2.5e (blue). 
For N = 10, the corresponding values are e w = 0.9e (green) and e w = 2.0e (blue). The substrate 
is 425<r long and width L y = A. 
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are shown at different times chosen so that the size of the precursor foot, when the foot 
spreads ahead of the droplet, is approximately the same size in each snapshot. In the weak 
wetting case (Fig. 2, left panels), the main droplet increases in size only slightly before 
saturating for A < 100a as shown in Fig. 3. For A = 30a, the bulk droplet saturates at a 
radius r& = 60a in both the wetting and non-wetting regimes. As A increases, the radius 
of the droplet in the non-wetting region decreases while the radius of the droplet in the 
wetting regime increases slightly. For A > 150(7, the main droplet spreads. Similar to the 
behavior seen in Fig. 1 for A = 200a, the A = 150a system in Fig. 2d (left) shows that 
the cross-section of material in the non-wetting region shrinks as the bulk droplet in the 
wetting region continually spreads. For A = 150a, the radius of the droplet in the wetting 
regime continues to grow at late time as shown in Fig. 3. Between A = 100a and 150a 
the gain in surface energy exceeds the cost of increasing the surface area and material is 
able to be drawn continually from the non-wetting region to the wetting region. Figure 2 
(right) presents N — 10 results for the strong wetting case. Compared to the weak wetting 
case where spreading of a precursor foot was nearly absent for A = 30a, Fig. 2a (right) 
clearly shows precursor spreading. A separate simulation of the strong wetting system with 
A = 20a showed a precursor foot spreads even for this very small feature size. It is also 
apparent, especially for smaller A, that the stronger interaction in the wetting region causes 
more spreading in the non-wetting region. The A = 150a case for strong wetting is not 
shown in Fig. 2; however, similar to the weak wetting case, the minimum width on which 
the bulk droplet spreads is between A = 100 and 150a. This indicates that the bulk wetting 
transition is related more to A than to the strength of interaction in the wetting region e w . 
For reference, A = oo is shown for the strong wetting case in Fig. 2d (right). 

N = 100 results for different A are shown in Fig. 4. At A = 60a, which is about 4 times 
the average end-to-end distance of the chain in the melt, the droplet spreads at early time to 
a radius = 80a before reaching its equilibrium size. This final size is comparable to what 
is observed for a homogenous surface with a mean surface interaction strength, e w ~ 1.75e. 
For A = 100a (Fig. 4b) and 200a (Fig. 1), the bulk radii in the wetting region saturate 
at lower values, m 70a. These data demonstrate that, for sufficiently small feature size, 
a mean interaction description is fairly accurate but becomes less so for increasing feature 
size. For A > 60a a precursor foot spreads ahead of the droplet while for A < 30a no foot is 
observed. The width A above which the main droplet spreads on the lyophilic strip is larger 
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FIG. 2: (color) Simulation snapshots for N = 10 polymer droplets; in all cases shown, the 
interaction strength of the wall in the non-wetting regions (green) is e w = 0.9e. For panels on the 
left, the interaction strength in the wetting regions (blue) is e w = 2.0e while on the right e w = 3.0e. 
In panels (a) - (c), snapshots are presented at t = 80, OOOr (left) and t = 30, OOOr (right) while in 
panel (d) t = 48, OOOr (left) and t = 15, OOOr (right). Pattern lengths are (a) A = 30a, (b) A = 60a, 
(c) A = lOOo", and (d) A = 150<r (left) or (d) A = oo (right). The substrate is 450<r long in each 
snapshot, while the width L y = A except for A = 30<r and A = oo, where L y = QOa. 

than we can access with present computational resources. For comparison purposes, Fig. 4c 
presents data for A = oo and it is obvious that the main droplet spreads for this limit of 
pattern size. 

The contact radius of the precursor foot increases as rj(t) = 2D eS t for all cases studied. 
The effective diffusion constant D e g normalized by its value for the homogeneous wetting 
case is shown in Fig. 5. A result observed in all cases is that D e g increases monotonically 
with A. Of note is an apparent cross-over in observed behavior with decreasing A for the 
N = 10 strong and weak wetting systems. Data imply that, at large pattern width, foot 
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FIG. 3: Radii of the bulk droplets for N = 10 in the center of the wetting regions for the weak 
wetting case e w = 0.9 and 2.0e for A = 30a (crosses), A = 60cr (diamonds), A = 100a (squares), 
A = 150o" (triangles), A = 200cx (circles), and A = oo (solid). 




FIG. 4: (color) Simulation snapshots for N = 100 polymer droplets with e w = l.Oe (green) and 
e w = 2.he (blue) for (a) A = 60a at t = 96, 000r, (b) A = 100a at t= 96, 000t, and (c) A = oo 
at t = 32, OOOr. The substrate is 425<r long in each snapshot, while the width L y = A except for 
A = oo, where L y = QOa. 

transport for the stronger wetting system is more affected by the presence of a lyophobic 
region than it is for the weaker wetting system. We hypothesize that this results from the 
very aggressive rate of transport for the foot on the stronger wetting substrate. The system 
must supply material to the foot from the bulk of the droplet. The rate of foot advancement 
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FIG. 5: Normalized effective diffusion coefficient for the precursor foot as a function of the surface 
pattern wavelength A for N = 10, e w = 0.9e and 2.0e (circles), N = 10, e w = 0.9e and 3.0e 
(squares), N = 100, e w = l.Oe and 2.5e (triangles). For Rq = 50a, the values for -D e ff,A=oo = 0.17 
(circles), 0.75 (squares), and 0.33a 2 jr (triangles). 

is large enough in the strong wetting case that it becomes limited by transport of material 

from the non- wetting region to the wetting region, and on to the foot. This is only true for 

large pattern width A > 60a where the foot is comprised of a more significant amount of 

polymer material. Note that foot transport along the lyophilic strip is observed for A < 30a 

in the strong wetting system whereas it is not for the weaker interaction. In our previous 

work studying wetting on homogeneous substrates, we found that D c s increased with 

i ii 

the size of the initial droplet as R . While we have not varied the size of the initial droplet 
in this study, we expect a similar dependence of D e s on the initial size of the droplet. 



IV. CONCLUSION 



Chemically patterning a surface offers a unique opportunity to control the spreading 
of liquid nanodroplets. Here we have presented simulation results which demonstrate the 
interplay between the pattern width A, the length of the polymer chain N, and the interaction 
strength of the lyophilic strips in determining spreading behavior. We find that there exists 
a minimum A, which is a function of N, below which spreading of the main droplet does 
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not occur. For N = 10, a bulk spreading transition with increasing A was observed, while 
for N = 100, the transition occurs for A larger than is presently accessible. In a number of 
cases, a precursor foot is observed to spread even when the main droplet does not spread. 
The observed effective diffusion constant of this precursor foot decreases monotonically with 
A. Future work will reveal the dependence of observed kinetics on the initial size of the 
droplet and specify the relation between N, A, and wetting transition. 
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